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RADIATIVE TRANSFER IN STELLAR ATMOSPHERES 

Jiff Kubat^ 



Abstract. This review presents basic equations for the solution of the 
NLTE radiative transfer problem for trace elements and methods for 

ms^ , its solution are summarized. The importance of frequency coupling in 

L^ ' radiative transfer in stellar atmospheres is emphasized. 

^' 1 Introduction 

o 

H ' A stellar atmosphere is the transition between the dense optically thick stellar 

^ , core and the transparent interstellar medium. It is the boundary layer of a star 

and it is the only part of the star we can see directly. Knowledge of the stellar 

atmosphere is a gate to the investigation of the stellar interior, and, consequently, 

to stellar evolution processes. The radiation passing through it serves not only as 

an information tool describing what is happening there, but it also significantly 

modifies the stellar atmosphere itself. This makes the task of determining the 

^4^ , emergent radiation from the star highly nonlinear and therefore also complicated. 

Cn| ' Since the stellar atmosphere is a transition from the optically thick stellar 

\f^ , interior to optically thin interstellar medium, the corresponding radiative transfer 

(^ • equation for the whole atmosphere cannot be subject to global simplifications. For 

^^ , the optically thick part of the atmosphere the diffusion approximation is valid. 

'""' ' Radiation propagates forward slowly resembling the process of diffusion. This 

^ , approximation greatly simplifies the radiation transfer equation. However, the 

J^ ■ diffusion approximation is unusable for the optically thin part of the atmosphere. 

I ^S , The situation there is again simple from another point of view. Since the medium 

C^ ' is optically thin, very little interaction with matter occurs there, which often allows 

us to neglect absorption and to take into account only emission. Between these two 

relatively simple regions there is a transitional region, where neither diffusion nor 

optically thin approximations are valid. We have to solve the full set of radiative 

transfer equations together with many constraint equations, which introduce the 

dependence of the absorption and emission coefficients on radiation and which 

make the problem to be solved highly nonlinear. 
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Radiation carries information about the medium where it was created. In stel- 
lar atmospheres, it has an additional role. Since the photon mean free path is much 
larger than the particle mean free path, radiation in stellar atmospheres influences 
the matter far from the place where it is formed, sometimes very significantly. It 
changes the population numbers of particular atomic energy levels, sometimes very 
far from the equilibrium values. These changes then influence spectral lines, which 
are observed by us. Radiation also transfers energy causing heating (or cooling) 
of specific parts of the stellar atmosphere. In some cases radiation has such an 
enormous effect on the atmospheric matter that it becomes the basic reason for a 
stellar wind, which is then called as 'radiatively driven'. All these crucial effects 
of radiation influence the transfer of radiation, which becomes strongly nonlinear. 

2 Model stellar atmosphere 

Construction of a stellar atmosphere model is a standard task of stellar atmo- 
sphere physics. It may be considered as a task to determine space distribution 
of basic macroscopic physical quantities (r is the position vector and v is the fre- 
quency), temperature T{r), electron number density ne{r), density /o(r), velocity 
v(r), radiation field J{i>,r), population numbers ni{r) of the level i, and others 
by solving equations for energy equilibrium (which determines T), radiative trans- 
fer (J), statistical equilibrium (n^), state equation (rig), continuity equation (p), 
equation of motion (v), and possibly also some other. For the simpler case of a 
static atmosphere in radiative equilibrium, equations of continuity and motion are 
replaced by the hydrostatic equation (which determines p) and the energy equa- 
tion simplifies to the equation of radiative equilibrium. Even in the static case 
the system of equations is quite complicated and further simplifications are often 
used. Once the atmospheric structure is known, detailed radiation field specific 
intensity I{i', n) for each direction n can be calculated by a simple formal solution 
of the radiative transfer equation. 

The final goal of stellar atmosphere modelling is straightforward. We want to 
compare the theoretical emergent radiation with that detected by observational 
instruments. In principle, it is possible to determine a detailed emergent radiation 
energy distribution from scratch (together with a model atmosphere). Due to the 
need to know many details of the emergent spectrum, both the number of vari- 
ables which need to be determined and consequently the number of equations to be 
solved become enormous. To make the problem tractable, this task is usually not 
solved at once, but it is divided into particular steps. First, the model atmosphere 
(preferrably NLTE) is calculated by taking into account all substantial contribu- 
tion of both physical processes and atomic data. Minor chemical elements and 
weak lines need not be taken into account in a great detail providing they do not 
influence the atmosphere structure significantly. Then the emergent radiation is 
calculated for a given model atmospheric structure taking into account all details 
of the chemical composition and line profile structure including the weakest lines. 

Alternatively, it is possible to insert an intermediate step after the model at- 
mosphere calculation, namely we may determine the NLTE level populations of 
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some ions not significantly influencing the atmospheric structure (trace elements) 
more accurately, and to solve the coupled set of radiative transfer and statistical 
equilibrium equations for them. For this restricted case with a given model at- 
mosphere (i.e. T and p are fixed) we do not solve the radiative and hydrostatic 
equilibrium equations and we solve only the radiative transfer and statistical equi- 
librium equations. This is the case of the NLTE line formation of trace elements 
in stellar atmospheres, which is the subject of these proceedings. In this case, 
we have to account for full angle and frequency dependence of radiation. In the 
subsequent step of emergent radiation calculation, they are considered as given. 
Note that in all steps the radiative transfer equation has to be solved. 



3 Basic equations 

3.1 Radiative transfer equation 

The full radiative transfer e quation for the specific intensity of radiation /(r , n, v, t) 



in three dimensions reads (jMihalasl . Il978l Eq. 2-24) 
1 dl{r,n,v,t) 



c dt 



(n • V)/(r, n, v, t) ^ -qir, n, v, t) - x{r, n, v, t)I{r, n, v, t). (3.1) 



Here r](r,n,i',t) is the emissivity and ^('"j'^i'^iO i^ the absorption coefficient 
(opacity). In some mathematicians' language this equation becomes a "7D" equa- 
tion (3 spatial dimensions, 2 directional cosines, time, and frequency). However, in 
the following, we shall use the notation "xD" referring only to the spatial dimen- 
sions. To describe fully the stellar atmosphere with all its properties and possible 
features (imagine a picture of the closest star - the Sun) , it is necessary to use the 
full 3D radiative transfer equation. This is a huge task and it may hardly give any 
reasonable results, even using contemporary fast computers. 

However, in some situations it is possible to use various simplifying assump- 
tions. The basic common assumption is the assumption of stationarity, which 
means that all d/dt — ;> 0, thus, all quantities are time independent. This assump- 
tion is used very often, however, it does not mean that it is necessarily always 
valid. Observational data are always averages over a certain time interval, which 
is usually much longer than the typical time scale for rapid changes in the atmo- 
sphere. This is, of course, true for distant objects, like most stars (where we have 
to collect enough light), but it fails for the closest star, our Sun. Nonetheless, in 
the following, we shall assume that the stationarity assumption is justified. 




Fig. 1. Schematic picture of the plane-parallel stellar atmosphere. Specific intensity of 
radiation depends on a z-coordinate, on the angle cosine n = cosS, and on frequency v. 
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If the stellar atmosphere is very thin compared to the stellar radius, then 
its curvature may be neglected and close parallel rays stay almost nearly parallel 
throughout the atmosphere. Then we can assume that the atmosphere is an infinite 
thin horizontally homogeneous plane. Such model stellar atmosphere is called 
plane-parallel (Fig.[T]). This is a very common assumption, which simplifies the 
geometry of a problem, it allows to use a ID radiative transfer equation instead 
of the full 3D one. Then, the radiative transfer equation takes the form (in the 
following we will often use the common notation convention I{v) — >■ I^) 

^dl^ = vAz) - xAz)I^u{z) = X. [SAz) - V(^)] , (3.2) 

dz 

where fj, — cos 9 is the angle cosine of the ray and 5',^ = rj^/xu is the source func- 
tion. Note that the one-dimensionality of equation (|3.2p means that the physical 
quantities depend only on one coordinate (z). The radiation field fills the whole 
space in and full angle dependence is taken into account, of course, with inherent 
symmetries, which allow to use only one angle cosine. 

If the stellar atmosphere is not thin compared to the stellar radius, the plane- 
parallel assumption may become invalid (e.g. for limb radiation calculations). 
Then another very common simplifying assumption, which also leads to one- 
dimensional atmosphere, may be used. It is the assumption of spherical sym- 
metry. This assumption is useful for the treatment of extended atmospheres. The 
radiative transfer equation p.ip then takes the form 

9V(r) ^ l-^^9V(r) 
dr r df-i 

Here r is the radial coordinate and ^ = cos 9 is again the angle cosine of the ray, 
as in the plane-parallel case. 

Alternatively, i f we introduce the so- called variable Eddington factor f^ = 
i /_^ IJ^Ifiu dfJ'/Ji' (JAuer &: Mihalasl . ll970l) . we may write this equation as a 2nd- 



l^^^^^ + / Z.' = V^r) - x.{r)I,Ar). (3.3) 



order differential equation for the first moment of the intensity J^, — 5 /^^ /^i/ d/i 
(the mean intensity), 

^^^[^4#^ = ^.(.)-^.(.), (3.4) 

where dr^ = — Xi/ dz is the optical depth. Opacity (xu) and emissivity (rj^) in the 
radiative transfer equation depend not only on local temperature and density, but 
also on radiation, which may propagate much longer distances than the mean free 
path of particles. This makes the problem non-local and non-linear. 

Equation (13. 4p has to be supplemented by boundary conditions, which for the 
case of the stellar atmosphere may be 

— — gy{z) ~ H:j for the upper boundary and (3.5a) 



dr^ 

d[^(z)J.(z) 
dr^ 



= H^ + g„{z) for the lower boundary, (3.5b) 
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where g^ = /g A* (/^^ ~ -^/Ti/) (^fJ-fJi^^ and H^ , H~ are incident fluxes at corre- 
sponding boundaries. It is common to assume zero incident radiation {H~ — 0) 
and diffusion approximation at the lower boundary (if+ — [l/3][di?i//dr,y]). 

Opacity and emissivity (and then the sou rce function 5',^) in stellar atmospheres 
can be calculated using the expressions (see iMihalad . 1 19781 Eqs. 7-1 and 7-2) 
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(3.7) 



The first term on the right hand sides of both equations p.6p and (|3.7p corre- 
sponds to line opacity (bound-bound transitions), the second one to bound- free 
transitions, the third one to free-free transitions, and the last one (only in l3.6|) to 
electron scattering. In these equations, n^, n;, and Uk are level populations (n* is 
the equilibrium value defined with respect to the ground level of the next higher 
ion), gi and gi are statistical weights, a denotes the cross section for the corre- 
sponding transition. Tig is the electron number density, and ae is the Thomson 
scattering cross section. Note that all level population may be heavily influenced 
by the radiation field. 

For the case where we can use the equilibrium values of population numbers 
(the so-called 'local thermodynamic equilibrium - LTE'), quantities n^ can be 
replaced by their equilibrium values n*, which depend only on local values of 
temperature T and electron density rig through the Saha ionization and Boltzmann 
excitation equilibrium laws. This greatly simplifies the problem, because both 
opacity and emissivity then depend only on local temperature and density, and 
are independent of radiation. Then the radiative transfer equation (|3.4p is linear 
in J,,. 



3.2 Equations of statistical equilibrium 

However, such a simplified description is not generally valid. In stellar atmo- 
spheres, the conditions may be far from thermodynamic equilibrium and Saha- 
Boltzmann equations for level population numbers cannot be used. In this case, 
level populations are determined using the equations of statistical equilibrium. 



Y, [ni {Ru + Cu)- n, {Ru + C,i)] = 0, 



i = 1,...NL 



(3.8) 
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This is the set of NL equations for each explicitly considered level. Since these 
equations are linearly dependent, an additional equation has to be used. Usually 
the charge or particle conservation equations are used. In the case of trace elements 
the equation which defines the abundance of the element with respect to some 
reference element (usually hydrogen, but not necessarily) may be used. 

The quantities Ru and Cu in the equation p.Sp are radiative and coUisional 
rates, respectively, for transitions between energy levels i and /. CoUisional rates, 

nidi = niUequiT), (3.9) 

where qu is a function of temperature only, do not depend on the radiation field 
(note that n*Cu = n*C'ii). On the other hand, both upward (Ru) and downward 
(Ru) radiative rates, 

n,Rii = n,47r / ^^^^ J^ dv (3.10a) 

niRu = ni^Air f ^|^ ('^ + J,j diy, (3.10b) 

depend explicitly on the mean radiation intens ity J^. Detailed expressions fo r 
radiative rates can be found in the next chapter ( Kubatl . l2010llthese proceedings! ) . 



3.2.1 Two-level atom 

The basic principles are usually best shown using very simple, so called textbook 
cases. The concept of a tw o-level atom fa lls within them. Let us assume an atom 
having only two levels (see lMihalasl . ll978 ). In this case the equations of statistical 



equilibrium p.Sp reduce to 

"2 {R21 + C21) - ni (i?i2 + C12) = (3.11) 

For this case it may be derived that the source function 

S, = {l-e) I ^,J,dv + eB,, (3.12) 

where ipi, is the line profile function, Bi, is the Planck function, and e = e'/(l+e')j 
and e' = C21 [1 — exp (— /ij^/fcT)] /^2i- Equation (|3.12l) clearly shows that the 
physical processes can be separated into the scattering ones, which depend on the 
radiation field J^, and to the thermal ones, which depend on the local Planck 
function. The lower the value of e, the less the LTE assumption is acceptable. 

3.3 Final radiative transfer equation witii a constraint of statistical equilibrium 



If we explicitly emphasize the quantity dependences in the equation (|3.2|) , it then 
takes the form 



dr.? 



J,(z)-S',(z,J,0 (3.13) 
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which shows its nonlinearity in the specific radiation intensity J. Using v' in- 
stead of V indicates that the radiative transfer in a particular frequency generally 
depends on the radiation field at all other frequencies. 

Equation (I3.13|) with the boundary conditions (I3.5P has to be solved using opac- 
ity p.6|) and emissivity (j3.7p . together with the statistical equilibrium equations 
(|3.8I) using both radiative p.lO|) and collisional (|3.9I) rates. These equations to- 
gether describe the basic NLTE radiative transfer problem which has to be solved 
in stellar atmospheres. 

4 Frequency coupling of radiation 

Let us turn our attention to another aspect of the radiative transfer problem, 
namely radiation coupling. One of the most striking features of radiation is the 
coupling of distant places, which would otherwise remain uncoupled. This coupling 
may occur over very long distances within stellar atmospheres. This long-distance 
coupling is described by the radiative transfer equation (I3.11l3.21l3.3p . Besides long 
distance coupling, scattering causes coupling between different radiation propaga- 
tion directions. 

An important aspect of radiation is the frequency coupling, which is crucial in 
stellar atmospheres. Although we succeeded to reduce the geometry complexity 
of the problem significantly, we have to retain full frequency dependence of the 
equations (13. 2[) and p.3p . since the radiative transfer in individual frequencies 
influences transfer in many other frequencies. Each of the equations p.2p or (|3.3p 
represents a set of frequency coupled equations. 

The trivial case is that of no interaction of radiation with matter. Then the spe- 
cific intensity / remains constant and there is no frequency coupling. The simplest 
case is that of pure coherent scattering, either continuum (like electron scattering) 
or in lines. In this case, photons may change direction, but not frequency and the 
radiative transfer equation has the form 

M^=/p-^/_^'^PP'Vd/i', (4.1) 

where cr^^' is the scattering cross-section. This kind of radiative transfer prob- 
lem forms a significa. nt part of the Chandrasekhar's Radiative Transfer book 
( Chandrasekhaii Il950l) . In this case the radiative transfer at individual frequencies 



may be solved independently of other frequencies. 

4.1 Frequency coupling within spectral lines 

In a spectral line, frequencies of absorbed and emitted photon may differ within 
the line profile, which for the case of pure natural broadening (Lorentz profile) is 

r 



2 



^^ = ^"\,,2 > (4-2) 



<"-""' +1,5 
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where vq is the central line frequency and F is the damping parameter. For the 
case of complete or partial frequency redistribution there is a coupling between 
them via the redistribution process. 

Thanks to thermal motions the spectral lines are subject to Doppler broaden- 
ing. As a consequence, even in the case of coherent scattering (which takes place 
in the atomic frame), thermal motions may cause the difference between absorbing 
and emitting frequency. The broadened line profile for the case of atomic frame 
coherent scattering (Doppler profile) is 



where the Doppler half- width Ai^d = (t'o/c) y^2kT/m. For the case of a naturally 
broadened line (|4.2p the profile is (Voigt profile) 

9?, = - ^H{a,v), (4.4) 



where H (a, v) is the Voigt function (see lMihaiaall978l Eq. 9-34). Thermal broad 



ening is dominant in stellar atmospheres, which strongly couples frequencies within 
each line. In addition other additional broadening mechanisms, line Stark or col- 
lisional broadening, may introduce further frequency coupling. 

To further illustrate the effect of lines on frequency coupling, we use the ex- 
ample of a two-level atom with partial frequency redistribution. The radiative 
transfer equation for this case is 



d/^ hv 
^~dl ~ 4^ 



-niip^BuJi, + UuiIju {Aui + Budu) , (4.5) 



where the emission profile ^{v) = J r (z/', ly) f(v') di^' , r {v' , v) is the redistribution 
function, ni and Uu are the occupation numbers of the lower (Z) and upper (u) 
levels, and A^i, B^i and Bi^ are the Einstein coefficients. All line frequencies, 
which are coupled via Equation (14.51) . have to be solved together. 

4.2 Frequency coupling across lines 

A photon after being absorbed, may be then re-emitted in the same line (this 
option was already mentioned in the preceding Section), or it may be emitted 
in a different line (or continuum transition), or it may be destroyed by collisional 
transition. In the latter case the photon's energy is converted to heat. As a reverse 
process, photons may be emitted after collisional excitation, which causes cooling. 
Each absorption or emission changes the distribution of atomic excitation 
states, which causes the change in the opacity. The changes are refiected by the 
set of statistical equilibrium equations (|3.8p . However, the radiative rates (|3.10p 
depend on radiation field. Consequently, we have to solve the equations p. 21) and 
simultaneously. 
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This coupling, which causes absorption of radiation in one part of the spectrum 
and its subsequent emission in a different spectral part, peaks in the effect usually 
referred to as line-blanketing. Spectral lines are not distributed across the fre- 
quency spectrum evenly. There appears a huge amount of lines in the UV spectral 
region. In addition, the radiative flux in this region is large. As a consequence, 
the radiation is absorbed in UV and re-emitted in the visual and infrared regions, 
where the opacity is lower. It is caused mostly by metallic lines dominated by Fe 
and Ni. The effect of line blanketing hcs been studied b y a number of authors 
Ce.g. iKurucd. I1979I: ICarbonl. [l98 4: Drcizlcr & Werner, 1993[ lHillier fc MiUeiJ . flOflJ 



Hubenv fc Land . Il995t iLanz fc Hubenvl . l2003l l2007l and references therein) 



4.3 Frequency coupling in moving atmospiieres 

For expanding atmospheres {v(r) 7^ 0, dv/dr > 0) another type of frequency 
coupling occurs. Since the process of atomic absorption and emission (scattering) 
takes place in a rest frame of an atom (i.e. in a co- moving frame), when it is 
observed from the reference frame connected with a star (the observer's frame), 
Doppler shift changes the frequency of lines. The radiative transfer equation in the 
observer frame has the form (|3.2p . but the opacity (Xtn^iz)) and emissivity {ri^^{z)) 
coefficients are angle dependent. The angle dependence is caused by the Doppler 
shift v' = V {1 — n ■ v/c) = v [1 — ^v / c). Equivalently, we may write the radiative 
transfer equation in the co-moving frame (we change the notation /_*' — > /i, v' — > v) 



M 



dz 



yrv dv 
c dz 



dlf,i,{z) 
dv 



Viy{z) -Xu{z)Ii,u{z). 



The second term on the left hand side clearly shows the frequency coupling caused 
by a derivative dl^i,{z)/di'. If we can neglect this term, we obtain the static 
radiative transfer equation. On the other hand, if we can neglec t the first term. 



we obtain the so-called Sobolev approximation f see ICastoii 120041 ) . 



5 Discretization 

Since the solution of the combined equations of radiative transfer and statistical 
equilibrium is performed numerically, we have to switch from continuous inde- 
pendent variables to a discrete set of representative points. Instead of the depth 
variable z we use depth points indexed by d. Since the density in stellar atmo- 
spheres varies by several orders of magnitude, the advantageous choice is to space 
the depth points equidistantly in log to (m is the column mass depth, dTO. — —pdz) 
or logr, where r is some suitable optical depth, e.g. Rosseland mean optical depth 
or an optical depth at some representative wavelength. The discretization of fre- 
quency V ^f n has to be done to resolve all spectral lines and continuum edges, 
consequently, it cannot be equidistant. On the other hand, for plane-parallel at- 
mospheres it is advantageous to use Gaussian quadrature for angle integration, 
which then defines the rays along which the radiative transfer equation is solved. 
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Discretized radiative transfer equation Discretizing the radiative tra nsfer 
equation (13. 4p for the frequency point n we obtain (see also iFeautrieiJ . |l96J) 

(5.1a) 
(5.1b) 
(5.1c) 



adfd-iJd-i + (bdfd + l)Jd + Cdfd+i'Jd+i = Sd for d == 2, . . . , ND - 1, 
{bifi+gi + l)Ji+cif2J2=H- +Su and 
ond/nd-i^nd-i + (^nd/nd + 3nd + 1) Jnd — H + ^nd- 
In these equations, 



ad = 

Cd = 

bd^ 

ci == 



AND 



K- 



jIAr,, 



-ad - Cd 



A7 



A' 



d+i 



At 



A^ 



1 -1 



d+i 



bi 



-ci, 



and 6nd 



-AND, 



(5.2a) 

(5.2b) 
(5.2c) 
(5. 2d) 

(5.2e) 



where At^_j. = Td — Td-i- The above equations refer to the case of 2nd-order 
differences, which is the most stable way of discretizing the differential ra diative 
trans fer equation. We may also use either splines or Hermite differences (see lAuer . 
19761 ) , which are more accurate, but sometimes we run into troubles since they are 
less stable and oscillatory behaviour occasionally occurs. 



Discretized frequency integration The integral both across the line profile 
and over continuum frequencies is replaced by a quadrature sum 



Jiv) Av 



NF 
n=l 



(5.3) 



NF 



For lines, it is necessary to ensure that J 0(i^) dz^ — '^n=i^n4'n = 1- If this 
condition is not fulfilled, it is obligatory to renormalize w„ . Otherwise, an artificial 
unphysical source or sink of photons appears, which may dramatically change 
results. 

To ensure correct treatment of ionization edges, it is advantageous to use two 
close frequency points which differ only by machine accuracy, one before and one 
after the ionization edge (Ranch, private communication). This ensures exact 
treatment of all contributions to the particular ionization rate and also to other 
rates where the continuum edge is a standard continuum frequency point. 



Discretized angle integration The integral over the angles in the calculation 
of the intensity moments is also replace by quadrature sums. 



» NF 



(5.4) 
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Usua lly, 3 points are sufHcient if Gaussian angle quadrature (see IChandrasekhai , 
ll95Cl Section 22) is used. Also in this case the check of correct normalisation is 
crucial. 

Discretized equations of statistical equilibrium For all d = 1,...,ND, 
i = 1 , . . . , NL we may write 

Y. {ini)d [iRH)d + {CH)d] - K)d [{Ra)d + {Cu)d]} = 0, (5.5) 

where the radiative rates are given by 

NF , . 

{ni)d{Rri)d = (nj)rf47rV'w„ ' " Jrf„, (5.6a) 

^ hVn 

n— 1 

NF ( \ / Oh 'i \ 

{ni)d{Ru)d = (nOd-4^E ^"^T^ i^ + Jdn). (5.6b) 

6 Formal solution of the radiative transfer equation 

Formal solution of the RTE is the next step after calculation of the model atmo- 
sphere and the occupation numbers. It is any radiative transfer solution for given 
Xv and r]i, (given S^) and it is relatively relatively simple. In this step the equation 
p.2p is directly solved, or it is sometimes rewritten to the 2nd-order form. The 
latter possibility has become more popular in one-dimensional model atmospheres. 
The formal solution is crucial for the total accuracy of the whole problem solution. 

6.1 Feautrier solution 

The most widely used radiative transfer equation solution for static pla ne-parallel 
atmospheres is the solution of the second-order system introduced by iFeautrier 



(|l964l ). In this scheme we combine opposite directions along a ray represented 
by specific intensities I"*" and /~. Introducing "Feautrier variables", i.e. their 
symmetric and anti-symmetric means. 







'^fiu — 


2 


(I^ 


+VJ 


and 
















Vf_ti, = 


1 
2 


(i;. 


+ VJ ' 


we obtain the transfer 


equation 












dr^ 


1/ 


= u^ 


U ~^ ^U 



3.1a) 



(6.1b) 



(6.2) 
supplemented by appropriate boundary conditions. 
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6.2 Short characteristics 

The short characteristics method is a first-order method. If we solve the radiative 
transfer equation for a finite slab (see Fig. [2]) with an optical thickness Ti, along a 
ray with an angle 6 {^ = cos 9), then by integration of Equation (|3.2|) between 
and Tiy we obtain for the outward direction 



and for the inward direction 



S'(i)e-(*-^-'/^dt, 



3.3a) 



C.(T.)=/,7,(0)e-^''/^ + 



(-m) Jo 



S'(i)e-(^--*)/''di. (6.3b) 



The stellar atmosphere is then divided into shells (i.e. finite slabs) between par- 
ticular depth points and solved by parts throughout the atmosphere. 

While not being used for the solution of the one-dimensional radiative transfer 
equation as much as the Feautrier solution, the power of the s hort characteris- 
tics m ethod emerges for multidimensional radiative transfer (see iKunasz fc Auer , 
19881) . 



7 A iteration 



We may formally write the formal solution of the equation p.4|) . which was de- 
scribed in th e preceding sec tion, as an operation of an operator A on the source 
function fsee lMihalasl . ll978h . 

J, = A,^,. (7.1) 

Generally, the source function, S^, depends on the radiation field intensity J^, 
which makes the problem non-linear. The simplest way how to cope with the non- 
linearity is an iterative procedure. Let us assume some starting values for opacity 
and emissivity, say the LTE values based on LTE level populations n* (i.e. the 
source function S^ = B^). The opacities and emissivities can be easily calculated 
from (|3.6p and p. 71) . respectively, using LTE population values. Then instead of 
the radiative transfer equation (|3.13l) . equation p.4p can be solved (for given rj^ 
and x^, i.e. for given S^), which is in fact the formal solution mentioned in the 



M 





i 



T, 



A- 



Fig. 2. Schematic representation of the short characteristics solution. 
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preceding section. Then, the equations of statistical equilibrium p.8|) are solved 
with radiative rates calculated for a given radiation field. New level populations 
rii are determined, which are then used for calculation of new values of Xv and t/j^, 
and the process is iterated. 

The above iteration scheme may be written as 

j(")=A,5("-i), (7.2) 

where (n) and (n— 1) denote the iteration numbers. The occupation numbers are 
determined using the equation 






{n|"^ [i?h(4")) + Ci^ - Ti^"' [i?,K4"^) + C.i] } = 0. (7.3) 



The level populations nf'' are used to calculate the source function SiP" and we 
follow again with solution of (|7.2I) , now for J^ ■ The crucial problem of this 
iterative process is that, for the case of stellar atmospheres, it converges extremely 
slowly, it rather stabilises far from the true solution. This is caused by the fact 
that each iteration represents interaction over one mean free path. In the optically 
thick parts of the atmosphere (with a lot of scattering events) it means extremely 
slow propagation of infor mation thro ughout the stellar atmosphere. Illustrative 
examples can be found in lAueiJ (jl984l ) . 



Although the lambda iteration is unusable for stellar atmospheres, for optically 
thin media, lik e planetary nebulae or circumstellar envelopes, it may work (see 



Dickel fc Aueij . 11994. ). However, for stellar atmospheres which are optically thick. 



different approaches have to be found. 

8 Complete linearization 

The complete linearization method (the Newton- Raphson metho d) was introduced 
to the theory of stellar atmospheres by lAuer fc MihalasI (|1968[ ) for the case of a 
full solution of a NLTE model plane-parallel atmosphere, which means simulta- 
neous solution of equations of radiative transfer, statistical equilibrium, radiative 
equilibrium, and hydrostatic equilibrium in one spatial dimension. For the sim- 
pler task, when the equations of hydrostatic and radiative equilibrium are not 
solved and only common solution of the radiative transfer equation together with 
the statistica l equi l ibrium equations is performed, the complete linearization was 
described bv lAueiJ (|l973h . 



We have to solve the plane-parallel radiative transfer equation p.4p together 
with the set of equations of statistical equilibrium p. 81) for radiation intensities Ji, 
and level populations n^. Since we cannot resolve the radiation field for infinite 
number of frequencies, we have to replace the continuous index (v) by a discrete 
one (i), J^ ^- Ji, and solve the radiative transfer equation only for these frequency 
points. The variables Ji and n; may be formally written as a vector 

tp = (Ji,..., Jnf; "•!,■• -j^-Nl) , 
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which has the dimension NF + NL. The equations to be solved may be formally 
written as 

F{^P)=0 (8.1) 



If ■0Q is the current estimate of the solution, which does not satisfy the Eq. (|8.1 
exactly, then the correct solution can be obtained by adding a correction 5xp, 



xp = ipQ + 5'4)- 



(8.2) 



Inserting it into the equation (|8.1|) we obtain after expansion to the first order for 
the corrections 

n -1 



Slj} 



dF 



(V'o) 



i^(^o)- 



(8.3) 



Then the new values of the vector tp are calculated using Eq. ()8.2I 



8.1 Standard formulation 



After implementing the difference scheme for the Eq. p.l3p . the linearized radia- 
tive transfer equation may be symbolically written as 



iSJd- 



id^Jd + CdSJd+l — Ld, 



(8.4) 



where the vector 5Jd ~ {SJds, ■ ■ ■ , ^Jd,NF) for each d = 1, . . . , ND. A^, B^, and 
<Cd are matrices, Ai — 0, and Cnd = 0. The right hand side Ld contains the 
source function Sd- 

The dependence of the source function on the occupation numbers may be 
expressed with the help of the derivative 






1=1 



(Sni), 



The combined equation may be then written as 

Ad5Jn,d-l + ^dSJn.d + CdSJn.d+l + ^dSUd ^ M nA- 



(8.5) 



If we formally write the statistical equilibrium equations p.8|) as An = b, where A 
is the rate matrix, n — (ni, . . . , riNL), we may then, for the changes of populations 
numbers Sud, write 



dA db 
on on 



6nri 



dA 



db 



aj" dj 



5Jd 



Ad ■ nd, 



which may be schematically written as 

EdSna + ¥d6Jd = K^. 



(8.6) 



(8.7) 
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Equations (|8.5p and (|8.7p schematically describe the problem, which we solve. 
Using the linearized equations of statistical equilibrium, we may express the level 
population changes as 



Sni,, 



NF 

E 



dni 
dJn 



SUn 



where 



r. NL 



dbr M—^ dA 



dJ, 



dJn 



(8.8) 



(8.9) 



Then the combined radiative transfer + statistical equilibrium equations may be 
written as a set of tridiagonal matrices 



where SJd ~ (SJi, ■ ■ ■ , 5Jnf), d — 1, 



(8.10) 



ND and A'^ = 0, C^^ = 0. 



8.2 Formulation using net radiative brackets 

Auer fc HeaslevI ( 19761) reformulated the application of the complete linearization 



method with the help of net radiative brackets. For each transition t between 
levels i and I, we may introduce the quantity 



Zt — UiRii — niRii, 



.11) 



which describes the net radiative rate. Linearizing this expression for each depth 
point d we obtain 



{SZt)^ = [ni)d{5Rii)d - {ni)d{5Rii)d. 



.12) 



Using Sui — ^j {drii/dZt) SZt we may eliminate the population numbers changes 
Srii and finally arrive at a system of equations 



-«'l|-^'li)--^'-«'gJ|--^ 



which is solved for 5[Zt)d- 



^^E^^^- (8-13) 



9 Accelerated lambda iteration 

Although complete linearization is a powerful method, it results in huge matri- 
ces having the dimension NF x NL. On the other hand , the simple A-i teration 
suffers from convergence problems. As an alternative, ICannonI ( 1973allbl ) sug- 
gested an iterative method base d on the operator split ting method. This method 
was later further developed by IScharmen ()l98ll Il984l) . Basic properties of this 
method, now called 'Accelerated Lambda Iteration' (ALI) method, were described 
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by Olson. Auer. fc Buchleij ( 19861 ). Useful reviews of ALI methods were published 
bv lHubenvl (Il99ll2003l) . 

We can introduce an approximate operator A* representing an approximate 
solution of the radiative transfer equation by adding and subtracting the term 
KlSv at the right hand side of (1711 . 

J,^KS, + {K-Kl)S, (9.1) 

In analogy with (j7.2p , we may introduce an iteration scheme 

4") = A^^t") + (A, - a:) 4"-i) = A:5(") + A4"-i) . (9.2) 

The basic difference with respect to the A-iteration is that the A*-operator is now 
acting on S{7 and not on S^ like in the case of ordinary A-iteration. The 

accurate solution of the radiative transfer equation A^S^ (the formal solution) is 
taken from the preceding, i.e. (n — 1) iteration step. In the limit of the converged 
solution the equation (19. 2p is exact. The term /S.Jv ' — [Ky — A*) S^ ' is the 
correction term, which describes, how the approximate solution (obtained with 
A*) for the current estimate of the source function S'("~^) differs from the exact 
one (obtained with Ai,). 

The normal procedure in multilevel radiative transfer calculations with ALI is 
to first determine the correct source function by common solution of the equations 
of statistical equilibrium with a simplified expression for 4" (j9.2p . and then to 
obtain the radiation intensity from the calculated source function using the formal 
solution. Since the source function is the ratio of emissivity and opacity, which 
both depend on population numbers, determination of the source function means 
in this case calculation of the population numbers of individual levels using the 
equations of statistical equilibrium. The expression (|9.2[) for the new estimate of 
the mean intensity is inserted to the expressions for the radiative rates p.lO|) . 



riiRa = nA-K 



niRii = ni—An 
91 



aii(y) 



Ksir^ + a4"~^^ 



dv 



aa(v) (2hv 
hv 



a:5(") + a4""^^ 



dv 



(9.3a) 
(9.3b) 



which are then used in the equations of statistical equilibrium (|7.3p . These equa- 
tions may be written as 



^{ni [Rii {ni,ni) + Cu] ~ riilRu {ni,ni) + C\i]} = 0, 



(9.4) 



i^i 



and they are solved for the population numbers n^. Finally, the formerly lin- 
ear set of the equations of statistical equilibrium became a nonlinear one. This 
nonlinearity is a tax for simplifying the problem. 

There are two basic possibilities how to solve this system of nonlinear equa- 
tions, namely the Newton- Raphson scheme (linearization), which was used, e.g., by 
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KubatI (|1994I ). or more popular approach of pre-conditioning (jRvbicki &: Hummeij . 
19911 ). The linearization scheme is very similar to that described in the Section [5] 
However, the savings achieved are very important, since we do not solve explicitly 
the problem in frequencies, number of which may reach several tens of thousands 
or more. 



9.1 Construction of the A* operator 

The crucial point is the construction of the operator A*. It should have two 
basic properties: First, it has to be simple enough to be calculated quickly, and, 
second, it should describe the most im portant interaction s of the radiation field 
with matter as accurately as possible. lOlson et al.l (|l986l ) showed that the best 
choice of such operator is simply the diagonal of the A operator. However, also 
other construction met hods of the ap proximate operator are being used. The 
tridiag qnal op e rator of lOlson fc Kunasz (1987) was used by Wern eij (|l989l ) and 
alsobv [Kubi!tl(|l994l) . 

The approximate lambda operator has to be numerically consistent wi th the 
formal solution of the radiative transfer equation, so the Olson fc Kunaszj ( 1987 ) 
operator is suitable if we use the short characteristic solution of the radiative 
transfer equation. On the other han d, for F e autrie r solution the approximate 
operator constructed after iRvbicki fc HummeiJ (|l99l[ ) has to be used. 



10 Summary 



In this introductory part we emphasized the importance of frequency coupling 
of radiation (in addition to the spatial and angle coupling) for the conditions of 
stellar atmospheres. This frequency coupling is ensured mainly by equations of 
statistical equilibrium, which connect distant parts of the frequency spectrum. 
We also briefly summarized underlying equations for the NLTE radiative transfer 
problem in stellar atmospheres with a focus on the problem of its solution for the 
case of trace elements. 

There are two basic methods of the solution of the RTE+ESE system, namely 



the complete linearizati on method, which is used by the Kiel code (see iKamnl . 
20ld Ithese proceedings), and the accelerated lambda iteration method used by 
the DETAIL code (see iButlen . I2010I Ithese proceedings! ) . 
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